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applied  to  the  fast  reaction  missile  defense  problem.  The 
results  are  compared  with  Proportional  Navigation.  The  states 
controlled  in  the  minimum  time  application  are  the  derivatives 
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I.  INTRODUCTION 

The  expectation  of  any  missile  system  is  to  hit  the 
target  with  high  probability;  in  other  words,  to  make  the 
terminal  miss  distance  between  the  missile  and  target  as  small 
as  possible.  For  thousands  of  years  projectile  weapons  have 
had  a  lack  of  control  after  launching.  This  is  detrimental  to 
hit  probability.  Guided  missiles  have  overcome  this  situation 
by  being  controlled  after  launching. 

Currently,  there  are  three  guidance  laws  for  tactical 
missiles.  These  laws  are  Pursuit  Guidance,  Line-of -Sight 
Guidance  and  Proportional  Guidance.  A  missile  can  employ  one 
or  more  of  these  guidance  laws.  Guidance  at  the  terminal  phase 
of  the  missile  flight  is  most  important. 

A  point  defense  system  may  be  required  to  hit  an 
attacking  missile  that  has  a  speed  advantage.  Instead  of 
following  an  expensive  method  such  as  designing  new  systems 
against  these  targets,  current  systems  may  be  upgraded. 

Present  work  applies  a  different  control  to  a  surface-to- 
air  missile.  In  Chapter  II,  a  minimum  time  control  problem  is 
defined  and  solved  for  a  third  order  regulator.  In  Chapter 
III,  the  three  different  missile  guidance  laws  are  described. 
Acceleration  effects  on  missile  velocity  is  also  explained.  In 


Chapter  IV,  a  missile  model  is  formed  and  the  missile/target 
engagement  scenario  is  simulated  by  using  Proportional 
Navigation  Guidance.  In  Chapter  V,  the  third  order  minimum 
time  solution  is  applied  to  the  missile/target  engagement 
simulation. 


n.  MINIMUM  TIME  CONTROL 

A.  OPTIMUM  CONTROL 

Control  system  design  is  a  trial  and  error  process 
utilizing  various  techniques  to  reach  a  desired  outcome. 
Classically,  the  performance  of  a  system  is  defined  in  terms 
of  time  and  frequency,  i.e.,  rise  time,  overshoot,  gain  and 
phase  margins.  However,  complex  multiple-input,  multiple- 
output  systems  require  a  more  elaborate  design  than  the 
classical  design.  For  example,  when  controlling  the  attitude 
of  a  satellite  to  minimize  fuel  consumption,  optimal  control 
theory  is  required  to  design  a  satisfactory  system, 
incorporating  the  complicities  of  multiple  input  and  output 
variables.  Optimal  control  theory  can  be  defined  as: 

"The  objective  of  optimal  control  theory  is  to  determine 
the  control  signals  that  will  cause  a  process  to  satisfy 
the  physical  constraints  and  at  the  same  time  minimize 
(or  maximize)  some  performance  criterion. " [Ref.  l:p.  1] 

B.  TIME  OPTIMAL  CONTROL 

1.    Problem  Definition 

In  minimum  time  problems,  the  objective  is  to  drive 
a  system  from  an  arbitrary  initial  state  to  a  desired  state  in 
minimum  time.  Mathematically,  the  problem  is  to  transfer  a 
system,  that  is, 


x(t)    =  Ax(t)    +  Bu(t)  (2.1) 

or, 

x(t)    =  f(x(t)  ,u(t)  ,  t)  (2.2) 

from  an  arbitrary  initial  state  to  a  desired  state  by 
minimizing 

J   =  jdt   =    tf  -    t0  (2.3) 

to 

Control  effort  u    is  constrained  by  a  maximum  value,  such  as 

\u\   ±  N  (2.4) 

The  final  state  may  be  any  point  in  the  state-space. 
In  this  study,  the  final  point  is  the  origin  of  the  state- 
space.  When  the  final  state  is  the  origin,  in  a  linear, 
stationary,  nth  order  system  as  defined  in  Equation  (1.1)  ,  the 
problem  is  referred  to  as  the  stationary,  linear  regulator, 
minimum-time  problem. 

2.    Problem  Solution 

In  this  section,  minimum  time  control  for  single 
input  systems  will  be  described.  From  Pontryaginf 1] ,  it  is 


known  that  the  performance  measurement  J  defined  at  Equation 
(2.3)  can  be  minimized  by  minimizing  the  Hamiltonian  which  is 
in  the  form  of 

H   =  1  +  pTAx  +  pTBu  (2.5) 

In  this  equation,  u  minimizes  the  Hamiltonian  by 
operating  at  its  maximum  value  with  the  opposite  sign  of  its 
coefficient  pTB.    Thus,  u   can  be  written  as 

u  =  -N*sign{pTB)  (2.6) 

An  optimal  time  control  can  be  solved  for  a  system, 
if  all  the  eigenvalues  of  A  have  non-positive  parts.  The 
control  u  is  called  the  bang-bang  control  which  occurs  at  its 
maximum  values  as  either  +N  or  -N.  The  number  of  switches 
between  these  values  is  at  most  n-1  times  for  systems  with 
real  roots  depending  on  the  initial  position  of  the  system  in 
state  space  where  n   is  the  order  of  the  system. 

For  an  existing  control,  there  are  two  trajectories 
possible  reaching  to  the  origin;  one  for  u=+N  and  one  for  u=- 
N.  These  are  zero  trajectory  curves  that  finally  carry  the 
system  to  the  origin.  These  trajectories  divide  the  state- 
space  into  two  parts. 


Some  possible  optimal  trajectories  for  a  second 
order  system  can  be  seen  in  Figure  2.1.  The  solid  line  curves 
are  the  zero  trajectory  curves.  The  state  of  the  system  above 
these  curves  follows  a  parabolic  path  parallel  to  the  -N  zero 
trajectory  curve,  with  the  control  effort  of  -N;  the  state 
below  these  curves  has  control  effort  +N  and  follows  a 
parabolic  path  parallel  to  the  +N  zero  trajectory.  These 
curves  are  illustrated  with  dashed  line  curves.  The  control  on 
these  curves  switches  once.  At  points  A  and  B  the  control 
switches  from  +N  to  -N  and  at  points  C  and  D  from  -N  to  +N. 
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If  these  trajectories  are  described  geometrically  in 
terms  of  state  variables,  the  control  function  can  be  defined 
as  the  switching  function.  This  control  function  decides  the 
sign  of  the  control  and  in  which  conditions  the  control  will 
change  its  sign  in  state-space. 

Since  the  final  state  is  the  origin,  but  the  initial 
state  is  any  place  in  state-space,  it  is  better  to  solve  the 
u  in  negative  time  with  zero  initial  conditions.  Zero  initial 
conditions  in  negative  time  define  the  origin,  end  point  of 
the  positive  time  problem.  In  this  method,  all  minimum  time 
trajectories  are  defined  in  a  generalized  equation  in  terms  of 
states.  This  process  is  explained  in  Reference  1. 

C.    TIME  OPTIMAL  CONTROL  EXAMPLE 

As  an  example,   switching  law  can  be  found  for  the 
following  third  order  system. 


U    1  l/s   X3   l/s     x2    l/s    Xi 

• ^- • »*• ♦ ^- • ^- ♦ 


This  system  is  similar  to  the  missile  model  which 
will  be  described  in  the  following  chapters,  so  the  switching 
law  which  will  be  solved  in  this  section  will  be  used  for  the 
minimum  time  application  of  the  missile  problem. 
The  state-space  of  the  system  is 


x. 


0    1 

0  " 

xi 

o" 

0    0 

1 

X2 

+ 

0 

0    0 

-a. 

X3. 

1 

u 


(2.7) 


The  transfer  function  of  the  system  is 


X(s)    _ 


U(s) 


s2  (s+a) 


(2.8) 


Eigenvalues  are  real  and  negative,  so  we  can  have  a 
solution.  Since  the  system  is  third  order,  there  are  at  most 
two  switches  in  the  control  depending  on  the  initial  position 
of  the  system  in  the  state-space. 

1.    Uncoupled  System 

The  switching  law  for  the  system  defines  hyper- 
surfaces  in  the  state-space  instead  of  curves  in  the  state- 
plane  for  a  second  order  system.  Transformation  of  the  system 
into  an  uncoupled  system  provides  an  easier  definition  for  the 
switching  law.  This  is  accomplished  by  taking  a  partial 
fraction  expansion  of  the  transfer  function  and  defining  each 
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state  with  respect  to  control  input  and  assigning  a  new  state 
vector  to  each  eigenvalue,  which  is 


x   =  Gy 


(2.9) 


so, 


Uncoupled  state-space  is  found  by  following 
Gy  =  AGy  +  Bu 

y  =  G^AGy   +  G^Bu  (2.10) 


Description  of  each  state  with  respect  to  control  is 


x1(s)    = 


u(s) 


S"  [S  +  OL 


1/g  _  1/oc2  +  I/a2 
s2  s  s  +  ct 


(2.11) 


Xo  (S) 


1/a    1/a 


s+a 


(2.12) 


x,  (s)    = 


s+a 


(2.13) 


From  these  equations,  G  and  G"1  matrices  are  found  as 


G  = 


_1 
a 

0 

0 


1 

1  ' 

a2 

a2 

1 

1 

a 

a 

0 

1 

and    G 


-l  _ 


a 

1 

o' 

0 

a 

1 

0 

0 

1. 

(2.14) 


The  uncoupled  state-space  system  is 


Yi 

0    10" 

Vi 

o' 

y* 

= 

0    0     0 

y2 

+ 

i 

A 

0    0    -a. 

y3. 

l 

(2.15) 


2.    Problem  Solution 

The  discretization  of  the  system  is 


<t>  =  ^((sI-A)'1)    = 


1  t   0  ' 

0  10 

0  0  e"Bt 


(2.16) 


A  =  f$Bdt   = 


li-2 


-±(l-e-at 
a 


(2.17) 


Using  these  discrete  time  matrices,  the  discrete  state-space 
equation  is 


>!<«" 

y2(t) 

= 

y3(t)_ 

1  t 

0  1 
0  0  e 


-at 


yx(o) 

y2(o) 
y3(o) 


2 


-  (l-e"at 
a 


u(0) 


(2.18) 
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The  positive  time  difference  equations  are 


yxU)    =  y^O)    +  y2(0)  t  +  ±u(0)  t2  (2.19) 


y2(t)   =  y2(0)  +  u(0)  t  (2.20) 


y3(t)  =  y3(0)e-at  +  iii°i  (i-e-at)        (2.21) 


In  negative  time  with  zero  initial  conditions  the 
difference  equations  become 

yx(t)    =   -|u(0)  t2  (2.22) 

y2(t)  =  -u(0)  t  (2.23) 

y  (t)  =  J£i£l(l-eat)  (2.24) 

a 

By  applying  u=+N  to  the  system,  we  have  Table  1  for 
one  second  intervals.  For  u=-N,  we  have  Table  2.  Zero 
trajectory  curves  are  defined  by  Tables  1  and  2. 

As  a  next  step,  the  curves  before  final  switching 
are  to  be  defined.  This  can  be  accomplished  by  equating 
positive  time  equations  to  the  negative  time  zero  trajectory 
curve  equations. 
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TABLE  1 
ZERO  TRAJECTORY  SET  (U=+N) 


u=+N 

t=i 

t=2 

t=3 

Yi(t) 

In 

2 

2N 

1" 

y2(t) 

-N 

-2N 

-2N 

y3(t) 

^(l-e«) 
a 

^(l-e2«) 
a 

^(l-e3«) 
a 

TABLE  2 
ZERO  TRAJECTORY  SET  (U=-N) 


u=-N 

t=i 

t=2 

t=3 

Yi(t) 

-In 

2 

2tf 

-2» 

2 

y2(t) 

tf 

2AT 

3N 

y3(t) 

-^<l-e«> 
a 

a 

a 
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At  the  time  control  switches  from  u=-N  to  u=+N, 
negative  time  equations  with  control  effort  u=+N  are  equal  to 
positive  time  equations  with  control  effort  u=-N.  These  sets 
of  equations  are 

yi(t)  =lN  =  yi(0)    +y2(0)  -  ±Nt2  (2.25) 

y2(t)    =-JV  =  y2(0)  -  Nt  (2.26) 

y3(t)  =  -(l-e«)  =  y3(0)e-at-^(l-e-at)      (2.27) 

At  time  t  the  control  switches  from  -N  to  +N.  From  Equation 
(2.26)  we  can  get 

t-l*^SL  (2.28) 

N 

By  putting  t  in  Equations  (2.27)  and  (2.28)  we  have 

*  =  yi(0)  *  Slj°L  (2.29) 

0  =  e-e"aZ^i(y3(0)+^\  -  M  +  i^e«       (2-30) 
\      a/    a    a 

By  following  the  same  steps,  we  can  solve  the  set  of 
equations  for  the  control  switches  from  +N  to  -N  and  obtain 
the  family  of  curves  shown  in  Table  3.  A  common  solution  from 
the  Table  3  gives  the  switching  law  for  the  uncoupled  system. 
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TABLE  3 
SOLUTION  SET  FOR  THE  UNCOUPLED  SYSTEM 


U=+N 

t=l 

-.*..).  is' 

-»y,«» 

0  -  e-e      *     (y,(0)    +  ■£)  -  M  ♦  i?e« 

\                  a  /         a         a 

t=2 

«■  *<°>-4£- 

-*ya(0) 
0  «  «-»■•      *     (y,(0)   ♦  -2)  -  M  ♦  i^e»« 

\                                  ft  /                ft                 tt 

t=3 

».*»..^ 

-«y,(o) 
0  =  es«e      *     (y,(0)    ♦  -J)  -  -^  ♦  -£eJ« 

u=-N 

t=l 

■*-*<« -*r 

0  ..-.     '»    (y,(0)   -*)*M-2?.. 

\  J               a/         a          a 

t-2 

.«.*,.>  -2^ 

•y,«H 

0  -  ea«e     *    ly.(O)   -  H)  +  M  -  He" 
\  3              a/         a         a 

t=3 

— *«•>-*£' 

0  .  .-.    '»    (y,(0)  -  »)  .  M  -  £e" 

\  J               a/         a         a 
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From  Table  3,  it  is  defined  ; 


•  J     ,m    y2(0)ly2(0) 
w  =  sigrly^Q)    +  -2 J^ 


(2.31) 


f  =  yx  ( o )  +  w 


y2(Q) 

2N 


(2.32) 


z  = 


N 


f 

AH 


(2.33) 


These  common  equations  define  the  family  of  curves  which  is 


"^^"^'^Trlf1?) 


(2.34) 


Equations  (2.31-34)  give  the  switching  law.  The  switching 
function  for  the  uncoupled  system  is 


u  =  -Nsign 


way,  (0) 

"— (y3<o)  ♦  iSW 


Y2  W) 


z2N+e^N 


a 


(2.35) 


This  switching  equation  is  for  the  uncoupled  system. 
For  the  normal  state  system  we  need  to  define  the  uncoupled 
states  in  terms  of  normal  states.  This  can  be  done  by  the 
following  equation 


y  =  G~xx 


(2.36) 
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So  the  uncoupled  states  are  expressed  as 

y1   =  ax1   +  x2 

y2   =  ax2   +  x3  (2.37) 

y3   =  *3 

From  Equations  (2.36)  and  (2.37),  the  switching  law  in  terms 
of  normal  states  can  be  written  as 

ii  =  -N*sigiL"m ■ (X3(o)+M)+v|^+^)j    (2.38) 

Where  z,w,f  are  also  in  terms  of  normal  states. 
3.    Simulation  of  The  Third  Order  System 

The  third  order  regulator  is  simulated  with  the 
switching  law  found  in  the  previous  section.  It  is  seen  that 
the  control  changes  its  sign  at  two  times  and  drives  the 
states  of  the  system  to  the  origin.  The  initial  values  for  the 
simulation  run  are  x^O^-0.5,  x2(0)=-0.5,  x3(0)=0.5.  Figures 
2.3-5  show  that  the  states  reach  to  zero.  The  control  effort 
is  in  Figure  2.6.  After  the  states  reach  to  the  origin,  the 
control  starts  chattering  or  limit  cycles.  Since  the 
simulation  is  in  discrete  time,  the  control  moves  back  and 
forth  between  plus  and  minus  zero  trajectories  around  the 
origin.  This  movement  causes  the  limit  cycles.  The  sampling 
rate  or  time  delay  in  control  effort  decides  the  magnitude  of 
the  limit  cycle.  Figures  2.6  and  2.7  show  the  second  order 
parabolic  relation  between  states  Xj  x2  and  x2  x3.  System  states 
follows  a  parabolic  path  in  the  second  order  state-plane. 
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Figure  2.2  First  State 
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Figure  2.3  Second  State 
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Figure  2.4  Third  State 


Figure  2.5  Control  Effort 
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Figure   2.6      State-Plane   for  x1  and  x2 


Figure   2.7      State-Plane   for  x2  and  x3 
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m.  MISSILE  GUIDANCE  AND  CONTROL 

A.  INTRODUCTION 

Guided  missiles  have  much  importance  over  other  weapons 
because  of  their  high  probability  of  hitting  their  target. 
This  high  probability  is  the  result  of  adjustments  in  missile 
movements  to  keep  the  missile  on  an  intercept  course  despite 
target  movements.  The  missile  acts  under  commands  from  the 
guidance  and  control  system.  These  systems  keep  the  missile  on 
an  intercept  course  with  the  target.  The  guidance  system 
applies  the  command  to  the  missile  body. 

B.  GUIDANCE 

An  unguided  weapon  may  have  an  excessive  miss  distance 
because  of  the  following  reasons: 

Incorrect  direction  of  launching; 
Perturbation  of  weapon  by  weather  or  wind; 
Unpredictable  movement  of  the  target  after  launching. 
A  good  method  to  reduce  miss  distance  and  thus  improve 
hit  probability,  is  to  use  a  closed  loop  system.  A  general 
form  of  a  closed  loop  system  is  shown  in  Figure  3.1  [Ref  2.  p. 
71] .  All  guidance  systems  are  particular  examples  of  the 
general  system  in  Figure  3.1. 


20 


Target. 
Data  I 


Guidance 
System 


Steering 
orders 


Missile 
Data 


Sensors  for 
Missile  Data 


Figure  3.1  Closed  Loop 

Missile  position  and  behavior  is  measured  by  an 
observation  instrumentation  which  may  be  actually  combined  in 
the  missile  itself  or  may  be  situated  on  a  remote  platform 
such  as  a  ship  or  an  aircraft.  This  missile  data  is  fed  into 
the  guidance  computer  with  target  data.  The  computer 
determines  the  maneuvers  of  the  missile  to  improve  its  chances 
of  hitting  the  target.  Steering  instructions,  such  as  desired 
lateral  acceleration  in  the  pitch  and  yaw  plane,  are  passed  by 
the  computer  to  the  control  system.  The  control  system  moves 
the  control  surfaces  or  determines  the  propulsive  thrust.  The 
resulting  motion  of  the  missile  is  measured  by  observation 
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instrumentation.  The  cycle  continues  with  new  observation 
data. 


C.    TYPES  OF  GUIDANCE  LAW 
1.    Pursuit  Guidance 

Pursuit  guidance  is  established  by  having  the 
missile  velocity  vector  directed  toward  the  target.  The 
missile  always  stays  along  the  line-of -sight  from  missile  to 
the  target.  This  guidance  is  effective  against  slow  moving  or 
oncoming  targets  which  have  slow  line-of-sight  rate;'  The 
pursuit  guidance  is  illustrated  in  Figure  3.2. 


Target  position  and  velocity  vector 
12      3      4 


^BBBBHBHI^^^I 


iigrimfr  iii 


Missile  position 
and  velocity  vector 


Figure  3.2   Pursuit  Guidance  Trajectory 
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2.  Line-Of -Sight  Guidance 

The  objective  of  a  Line-of -Sight  guidance  system  is 
to  constrain  the  missile  to  lie  as  nearly  as  possible  on  the 
line  joining  the  observation  point  and  the  target.  This  may  be 
done  by  two  forms  which  are  called  Command  Guidance  to  Line- 
of  -Sight  and  Beam  Riding.  In  beam  riding  guidance,  the  target 
tracker  maintains  its  radar  beam  on  the  line-of-sight  and  the 
missile  tries  to  stay  on  it.  In  command  guidance  to  line-of- 
sight,  the  missile  receives  commands  to  stay  on  the  beam.  In 
Figure  3.3,  beam  riding  line-of-sight  is  illustrated  [Ref.  3]. 

3.  Proportional  Navigation  Guidance 

In  this  guidance  law,  the  missile  moves  in  such  a 
way  that  its  rate  of  turn  is  proportional  to  the  rate  of  turn 
of  line-of-sight  from  the  missile  to  the  target.  The 
navigation  ratio  is  a  fixed  or  variable  ratio  between  the 
missile  rate  of  turn  and  the  rate  of  turn  of  the  line-of- 
sight.  This  guidance  system  demands  an  acceleration  command 
perpendicular  to  the  line-of-sight.  Since  proportional 
navigation  anticipates  the  future  position  of  the  target,  it 
is  effective  against  maneuvering  targets.  The  proportional 
navigation  guidance  is  illustrated  in  Figure  3.4  [Ref 3]. 

D.    CONTROL  SYSTEM 

The  control  system  has  two  important  functions;  the  first 
one  is  to  provide  a  stable  flight  in  all  phases  of  flight  and 
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Figure  3.3  Line  of  Sight  Guidance  Trajectory 

in  every  probable  disturbance.  For  example,  in  the  presence  of 
wind,  the  necessary  balancing  commands  to  the  control  surfaces 
come  from  the  control  system.  The  second  function  is  to  follow 
the  acceleration  commands  coming  from  the  guidance  system  to 
keep  the  missile  on  an  intercept  course  with  the  target. 

The  acceleration  commands  are  provided  by  the  guidance 
system  and  applied  to  the  missile  by  the  autopilot  in  the 
control  system.  In  Figure  3.5a  [Ref.  4]  an  acceleration 
command  in  j  direction  in  missile  coordinate  system  is  seen 
for  a  missile  with  a  velocity  in  i  direction.  The  angle  7  is 
measured  from  the  reference  axis  shown  in  the  figure.  Figure 
3.5b  [Ref.  4]  is  the  top  view  of  ij  plane  and  shows  the  effect 
of  centripetal  acceleration  Aj  in  the  direction  of  missile 
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Figure  3.4  Proportional  Navigation 

velocity  v.  The  change  of  velocity  vector  along  j  axis  can  be 
seen  in  Figure  3.5c.  From  Figure  3.5c  [Ref.  4],  the  following 
equation  may  be  written  by  assuming  the  change  of  direction 
happened  in  At  seconds. 

Aj   =  V-|l  =  V?  (3.1) 

In  Equation  (3.1)  ,  commanded  acceleration  AjC  is  shown  in  terms 
of  tangentional  velocity  v   and  angular  rates  7 . 
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Figure  3.5  Central  Acceleration  Effect 
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IV.  MISSILE  MODEL  AND  SIMULATION 

A.  ASSUMPTIONS 

For  the  missile/target  engagement  scenario  it  is 
necessary  to  have  a  missile  model.  In  the  simulation,  the 
missile  will  be  considered  a  point  of  mass.  The  simplified 
missile  motion  equations  will  be  defined  under  the  following 

assumptions. 

The  seeker  head  angle  rate  is  an  estimate  of  the  line- 
of -sight  rate. 

The  missile  thrust  cancels  drag. 

The  reference  for  the  flight  path  angles  is  the  x   and 
y  plane. 

B.  TWO  DIMENSIONAL  MISSILE  TARGET  GEOMETRY 

Before  looking  to  the  missile  model,  it  is  advantageous 
to  understand  the  equations  of  motion  in  two-dimensional 
geometry.  In  Figure  4.1,  a  missile  and  target  geometry  is 
shown  in  x  and  y  directions.  VM  and  VT  are  the  velocity 
vectors  of  missile  and  target,  respectively.  7M  and  yT  are  the 
velocity  angles  from  the  reference  line.  The  missile  flight 
path  angle  can  be  expressed  as 


yM  =  Arc  tan 
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Figure  4.1  Intercept  Geometry  For  Two  Dimensions 


Where  V,^  is  the  velocity  component  in  the  x  direction  and  V,^ 
is  the  velocity  component  in  the  y  direction.  From  the  same 
geometry,  the  target  flight  path  can  be  expressed  as 


yT  =  Arc  tan 


(    V. 


TY 


TX 


(4.2) 


Where  V-^  and  V^  are  the  target  velocity  components.  The  range 
between  the  missile  and  the  target  is  defined  as 
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R  =  yj  (xT-xM)2  I  (yT-yM)2  (4-3) 

Where  xT,  yT,  xM/  yM  are  the  x   and  y  coordinates  of  the  target 
and  the  missile. 

a   is  the  angle  of  line-of-sight  from  the  missile  to  the 
target.  It  is  defined  as 


o  =  Arc  tan 


(4.4) 


Vt  ~  Ym 


\     XT    XM 


The  magnitude  of  velocity  vectors  can  be  defined  as 


yM  =  v'  ^  +  ^ 


(4.5) 


VT    =  yj     V^     +  V^y 


Line-of-sight  rate  (a)  is  a  necessary  factor  in 
proportional  navigation.  This  data  is  provided  by  the  seeker 
head.  Analytic  expression  for  the  line-of-sight  rate  is  : 

6  =  Ur~*„)  (Vjy-Vt/y)    -  JYt-Ym)  (vtx~vmJ  (4#6) 
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C.    MISSILE  MODEL 

1.  System  Signal  Flow  Graph 

The  total  system  signal  flow  graph  in  the  x 
direction  is  illustrated  in  Figure  4.2  [Ref.  4].  Since 
proportional  navigation  will  be  simulated,  the  guidance  system 
has  a  seeker  head  which  provides  line-of -sight  rate,  (a) .  The 
input  for  the  autopilot  is  the  multiplication  of  the 
navigation  ratio  and  a.  The  output  signals  of  the  autopilot 
are  required  signals  for  the  missile  dynamics. 

2.  Guidance  System 

In  signal  flow  graph,  the  guidance  system  includes 
the  seeker  head.  A  seeker  head  is  a  homing  head  mounted  on  an 
airborne  missile.  The  purpose  of  the  seeker  is  the  tracking  of 
the  target.  This  is  done  by  sensing  radiation  or  reflected 
energy  from  the  target.  Tracking  the  target  by  seeker  shows 
the  angular  direction  to  the  target. 

The  signal  flow  graph  of  a  seeker  which  is  gimballed 
to  the  missile  is  shown  in  Figure  4.2.  This  seeker  points  at 
the  target  by  rotating.  Applied  torques  provide  rotation 
proportional  to  the  target  of f bore-sight .  This  can  be  written 

T=J|3  (4.7) 

Where  T  is  the  torque,  J  is  the  moment  of  inertia  of  the 
seeker  head  and  J9  is  angular  acceleration.  From  Figure  4.2  and 
Equation  (4.7)  the  resulting  equation  of  motion  is 
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Figure  4.2  System  Flow  Graph 
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0  =  1  =  -K2£   -  Kx$   +  KXQ  (4.8) 


Where  Kj  and  K2  determine  the  time  constants  of  the  seeker.  The 
transfer  function  of  the  seeker  head  from  input  line-of -sight 
angle  (a)  to  the  output  seeker  angle  (/9)  can  be  expressed  as 


P(s)  ,  * 


i 


o(s)        s2   +  iqs  +  jq 


(4.9) 


The  flow  graph  of  the  seeker  transfer  function  is  shown  in 
Figure  4.2. 

3.  Flight  Control  System 

Previously,  it  was  mentioned  that  the  control  system 
provides  a  stable  flight  to  the  missile,  and  passes  the 
guidance  system  commands  to  the  control  surfaces  by  keeping 
the  missile  stable.  These  stability  requirements  are  satisfied 
by  the  autopilot  in  the  control  system. 

The  control  system  causes  a  lag  in  the  missile 
system.  In  the  present  missile  simulation  model,  the  autopilot 
is  modeled  as  a  first  order  lag  with  a  time  constant,  1/a. 

4.  Missile  Dynamics 

The  missile  can  be  considered  a  point  mass  moving 
under  the  acceleration  commands  perpendicular  to  the  missile 
velocity  vector.  From  Figure  4.3,  components  of  acceleration 
can  be  expressed  as: 
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Figure  4.3  Missile  Acceleration  and  Components 

aMx   =  -<3„sinY 

aMy    =   <3«COSY 


(4.10) 


The  position  of  the  missile  in  the  x  and  y  plane  can 
be  found  by  integrating  the  acceleration  components  twice. 
Figure  4.4  shows  the  missile  dynamics  in  two  dimensions. 

D.    PROPORTIONAL  NAVIGATION  SIMULATION 

In  this  section,  the  missile/target  engagement  will  be 
presented  in  two  different  simulations  by  using  proportional 
simulation.  In  the  first  simulation  the  defensive  missile  has 
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Figure  4.4  Missile  Dynamics  in  Two  Dimensions 

the  speed  advantage;  in  the  second  simulation  the  target 

missile  has  the  speed  advantage.  The  following  assumptions  are 

made  for  both  simulations 

The  missile  is  limited  to  20  g's  acceleration. 

The  target  is  on  the  final  leg  of  its  flight  and  is  now 
on  a  non-maneuvering  trajectory,  so  it  has  no 
acceleration. 

The  x  and  y  plane  is  the  reference  system  for  the  flight 
paths  and  angles  of  the  target  and  the  missile. 

The  miss  distance  will  be  estimated  by  interpolating  the 
two  most  minimum  ranges  between  the  missile  and  the 
target. 

-  The  navigation  ratio  is  3.25.  Flight  constant  time 
constant  is  1/3  second. 

K]  and  K2  constants  in  seeker  head  are  100  and  20, 
respectively. 
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1.    Speed  Advantage  of  Missile 

In  the  first  scenario,  the  defensive  missile  has  the 
speed  advantage  over  the  incoming  target.  The  initial 
conditions  for  the  simulation  are  : 

vmx(°)  ~  2000  feet/sec 

xT(0)   =  10,000  feet 

yT(0)   =  1000  feet 

vtx(°)  =  -1000  feet/sec 
All  other  initial  conditions  are  zero.  As  shown  in  Figure  4.5, 
a  successful  intercept  occurs  and  the  range  goes  to  zero.  Miss 
distance  is  0.625  feet.  Figure  4.6  is  a  plot  for  the  missile 
and  the  target  flight  paths.  The  target  starts  its  motion  from 
10,000  feet  distance  and  moves  towards  the  y  axis.  The  missile 
has  an  intercept  path  with  the  target.  Figure  4.7  is  a  plot 
for  the  applied  missile  acceleration.  Missile  acceleration 
components  in  x  and  y  direction  are  shown  in  Figure  4.8.  At 
initial  phase,  maximum  acceleration  especially  in  y  direction, 
is  commanded.  The  line-of -sight  angle  is  in  Figure  4.9.  Line- 
of-sight  angle  rate  is  shown  in  Figure  4.10.  Initially  these 
parameters  increase,  then  decrease  until  just  prior  to 
intercept. 
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Figure  4.7  Applied  Acceleration 
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Figure  4.10  Line-of-Sight  Angle  Rate 
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2 .    Speed  Advantage  of  Target 

In  this  simulation,  the  incoming  target  has  the 
speed  advantage  over  the  defensive  missile.  Only  the  missile 
and  target  velocities  are  changed  from  the  previous 
simulation.  The  velocities  are: 

VnJO)  =  1000  feet/sec 

VTx(°)  =  -2000  feet/sec 
All  other  initial  conditions  are  the  same  as  in  the  previous 
simulation.  As  seen  in  Figure  4.11,  the  range  goes  towards 
zero  but  with  an  increasing  miss  distance.  The  miss  distance 
is  17.96  feet.  This  is  an  expected  result  because  of  the  speed 
advantage  of  target.  Figure  4.12  shows  the  missile  and  the 
target  flight  paths.  The  missile  still  has  the  capability  to 
follow  the  target  motion.  This  is  due  to  present  geometry 
advantage  for  the  missile  and  non-maneuvering  target  motion. 
Applied  acceleration  and  acceleration  components  are  shown  in 
Figure  4.13  and  Figure  4.14,  respectively.  The  acceleration 
especially  in  the  y  direction  increases  during  the  initial 
phase.  Figure  4.15  is  a  plot  of  the  line-of-sight  angle.  Line- 
of-sight  rate  is  shown  in  Figure  4.16.  At  initial  phase  line- 
of-sight  rate  increases  rapidly,  then  decreases  as  the  missile 
comes  to  a  constant  course. 
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Figure  4.12  Missile  and  Target  Trajectories 
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Figure  4.13  Applied  Acceleration 
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Figure  4.14  Acceleration  Components 
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Figure  4.15  Line-of-Sight  Angle 


Figure  4.16  Line-of-Sight  Rate 
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V.  MISSILE  APPLICATION  OF  THIRD  ORDER  CONTROLLER 

A.    MISSILE  AND  THIRD  ORDER  REGULATOR 

When  the  system  dynamics  in  Figure  4.2  are  examined,  the 
autopilot  and  the  missile  dynamics  together  form  a  third  order 
system.  This  third  order  system  is  basically  similar  to  the 
third  order  regulator  switching  function  found  in  Chapter  II. 
The  difference  between  the  two  systems  is  the  factor  of  Vmcos7 
and  VMsiny    between  7  and  the  acceleration  components. 

The  autopilot  and  the  missile  dynamics  together  can  be 
considered  a  third  order  regulator  that  drives  the  selected 
states  to  zero  in  minimum  time  as  illustrated  in  Figure  4.1. 
If  the  right  states  are  chosen  to  be  driven  to  zero  in  state- 
space,  the  intercept  occurs  between  the  missile  and  the 
target. 

The  system  signal  flow  graph  of  the  missile  for  minimum 

time  application  is  shown  in  Figure  5.1.  In  simulation,  the 

same  initial  conditions  used  in  Chapter  IV  for  the  speed 

advantage  of  the  target  simulation,  will  be  used.  The  states 

which  will  be  driven  to  zero  are  : 

The  vertical  distance  between  the  missile  and  target  and 
its  derivative; 

The  derivatives  of  line-of-sight  angle  (a) . 


43 


D) 

a 

G 

•H 

O 

XI 

•rH 

V 

■P 

4-> 

U 

•H 

rt 

£ 

p 

CO 

fa 

iH 

1 

O 

CO 

M 

CD 

4-> 

4-> 

C 

fC 

O 

4-) 

O  W 


Figure  5.1  System  Signal  Flow  Graph 
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B.    GUIDANCE  BY  VERTICAL  DISTANCE  BETWEEN  MISSILE  AND  TARGET 

In  the  first  approach,  the  states  to  be  driven  to  zero 
are  the  vertical  distance  between  the  missile  and  the  target 
and  the  first  and  second  derivatives.  The  purpose  is  to  force 
the  missile  to  rise  to  the  height  of  the  target  in  minimum 
time,  then  turn  on  the  target  for  head-on-collision. 

Since  the  VMcos-y  is  a  multiplier  in  y  channel  of  the  flow 
graph,  it  is  taken  as  a  scaling  factor  for  finding  the 
uncoupled  system.  The  y  channel  of  the  system  can  be 
considered  as  the  following  third  order  regulator: 


Vitcosy 

u 

1 

l/s 

*{^Y    l/s 

y 

i/s  y 

-a 

/**                        *3 

Xa 

Xi 

Figure  5.2  Missile  y  Channel  as  Third  Order  Regulator 


The  differential  equations  of  the  system  are 


y  =  kx   =  x2  (5.1) 


y  =  x2   =  x3  (5.2) 
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y  =  x2   =  V*4  =  V(-ax4   +  u) 
=  v(-a^+u) 


(5.3) 


=  -ox3  +  uV 


Where  V=VMcos-y.    So  the  state-space  of  the  system  for  y  channel 
is 


y 

*i 

0     1       0  ' 

*i 

0' 

y 

= 

*2 

= 

0    0     1 

x2 

+ 

0 

x 

k\ 

0   0    -a. 

X3. 

V 

u 


(5.4) 


For  the  uncoupled  state-space,  the  states  are  defined  as 


V/g    _    V/a2    +    V/a: 


s+a 


(5.5) 


V/a        V/a 


s+a 


(5.6) 


*3    = 


V 


s+a 


(5.7) 


G   and   G"1   matrices   are 


G  = 


V/a  -V/a2  V/a2 
0  V/a  -V/a 
0  0  V 


aiid      G 


■1  _ 


a/V  l/V     0 
0      a/V  1/V 
0        0      1/V 


(5.8) 
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The  uncoupled  state-space  is 


Vi 

0    10" 

Vi 

o' 

y2 

= 

0    0     0 

y2 

+ 

1 

y3. 

0   0    -a. 

y3; 

1. 

u 


(5.9) 


It  is  seen  that  the  uncoupled  system  is  the  same  for  the 
third  order  regulator  example  in  Chapter  II.  Thus,  the 
switching  law  is  the  same  as  in  Equation  (2.35)  for  the 
uncoupled  state-space.  But  the  transformation  from  uncoupled 
state-space  to  normal  state-space  by  Equation  (2.36)  will  be 
different.  The  uncoupled  states  in  terms  of  normal  states 
are 


Vi    = 

—  X,     +    —  X, 

V    1         V   2 

y2  = 

—  x,    +    —  X, 
V    2         V   3 

y3  = 

lx3   =y 

V   3         ' 

(5.10) 


It  is  seen  that  all  states  are  divided  by  VMcos-y. 
Therefore,  the  switching  law  from  Equation  (2.38)  is  obtained 
by  dividing  the  states  by  VMcos7. 

After  running  the  simulation,  the  missile  rises  to  the 
height  of  the  target  and  then  turns  onto  the  target  for  head- 
on-collision.  The  flight  paths  of  the  missile  and  the  target 
are  in  Figure  5.3.  The  range  according  to  time  is  shown  in 
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Figure  5.4.  The  miss  distance  is  0.07  feet.  The  vertical 
distance  between  the  missile  and  the  target  according  to  time 
is  shown  in  Figure  5.5.  It  is  seen  that  the  vertical  distance 
is  driven  to  zero  and  kept  at  that  value  until  intercept 
occurs.  The  acceleration  and  the  acceleration  components  are 
shown  in  Figures  5.6  and  5.7,  respectively.  Acceleration, 
especially  at  the  beginning  of  flight,  is  applied  at  its 
maximum.  Acceleration  in  y  direction  starts  with  a  high  value 
and  decreases  as  the  vertical  distance  decreases. 
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Figure  5.3  Missile  and  Target  Trajectories 
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Figure  5.4  Range  vs  Time 


Figure  5.5  Vertical  Distance  vs  Time 
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Figure  5.6  Applied  Acceleration 


Figure  5.7  Acceleration  Components 
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C.    GUIDANCE  BY  DERIVATIVES  OF  a 

The  second  group  of  states  which  is  chosen  are  the 
derivatives  of  line-of -sight  angle  (a ,  a ,  a) .  If  the  a  is  driven 
to  zero,  the  a  becomes  constant.  If  the  relative  velocity  is 
decreasing  and  the  line-of -sight  is  constant,  intercept  occurs 
between  the  missile  and  target.  This  is  the  theory  of 
Proportional  Navigation  Guidance  and  illustrated  in  Figure 
3.4. 

We  can  have  the  a  and  its  derivatives  analytically  from 
Figure  4.1.  In  an  actual  application,  these  values  can  be 
estimated  by  a  Kalman  filter  or  a  Luenberger  observer.  Line- 
of-sight  angle  is  given  by  Equation  (4.4).  Derivatives  of  a 
are 

6  =  XVy  ~  yV*  (5.11) 


g  =  xay  ~  y*x   _    2  {xVx   +  yVy)  (xVy   -  yVx)  (5.i2) 


...  _  ( V^y+xy-Vya^yx)  (R2)  -2  {xVx+yVy)  (xay-yax) 

i?4 


2  (vl+xa^Vy+yay)  (xVy-yVx) +2  (xay-yax)  (xVx+yVy)   (5.13) 


RA 

x'J'y'      v-^  ry  J  w  xl 


8  (xV+yVv)  2  (xV-yVx) 
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Where  x  =  xT-xM  y  =  yT-yM 

Vx  =  V^-V^  Vy  =  V^-V^ 

ax    =     aTx~aMx  ay    =     aTy~aMy 

After  the  simulation,  we  found  that  the  control  drives 
the  derivatives  of  a  to  zero  and  intercept  occurs.  Flight 
paths  of  the  missile  and  the  target  are  shown  in  Figure  5.9. 
The  range  is  plotted  in  Figure  5.10.  Miss  distance  is  0.045 
feet.  The  line-of-sight  angle  is  shown  in  Figure  5.11.  Figure 
5.12  shows  the  line-of-sight  angle  rate.  It  is  seen  that  the 
control  drives  the  line-of-sight  angle  rate  to  zero,  so  the 
line-of-sight  is  kept  constant  and  intercept  occurs. 
Acceleration  and  its  components  are  shown  in  Figures  5.13  and 
5.14. 
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Figure  5.10  Line-Of -Sight  Angle  vs  Time 
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Figure  5.12  Applied  Acceleration 
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VI.  CONCLUSIONS 

A  third  order  minimum  time  switching  law  was  developed. 
The  switching  law  that  was  found,  was  applied  to  the  third 
order  fast  reaction  missile  defense  problem.  The  Proportional 
Navigation  Guidance  is  applied  to  the  same  problem  for 
comparison. 

Minimum  time  application  was  established  using  two 
different  approaches.  One  was  with  the  vertical  distance  and 
derivatives,  the  other  was  the  first,  second  and  third  order 
derivatives  of  Line-of -Sight  angle.  Both  approaches  gave  more 
effective  results  than  Proportional  Navigation  Guidance. 

Since  bang  bang  control  always  uses  negative  or  positive 
maximum  control  effort,  another  control  logic  is  necessary  to 
shut  off  the  control  effort  when  the  desired  conditions  are 
met.  Otherwise,  the  system  starts  a  chattering  mode  or  limit 
cycles.  This  may  be  a  subject  for  future  projects. 
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APPENDIX 

1.  PROPNAV.M 

%  Missile/Target  engagement  simulation 

%  Proportional  Guidance 

%  Target  has  a  non-maneuvering  flight 

%  Program  terminates  when  the  range  starts  to  increase 

%  Miss  Distance  is  found  by  interpolating  the  minimum  range 


clear, 

tclg 

Nr=3.25; 

alp=3, 

• 
f 

maxac= 

=640; 

Am=[0 

1  0 

o; 

0 

0  0 

o; 

0 

0  0 

i; 

0 

0  0 

o]; 

Bm=[0 

0;l 

0 ;  0  0 ;  0  1  ]  ; 

Aa=[-alp] ; 

Ba=[alp]; 
As= [0     1 ; 

-100   -20]; 
Bs=[  0;100]; 


%  Navigation  Ratio 
%  1/a  =  Autopilot  time  constant 
%  Allowed  maximum  acceleration 
%  Missile  state  matrices 


%  Autopilot  states 


%  Seeker  states 


%  DISCRETIZE  THE  SYSTEM 

dt=0.01; 
[phim,delm]=c2d(Am,Bm,dt) ; 

[phit]=[phim] ; [delt]=[delm] ; 

[phia,dela]=c2d(Aa,Ba,dt) ; 

[phis/dels]=c2d (As,Bs,dt) ; 

%  INITIALIZE  THE  STATES 

xm( : , 1)=[    0       %  xmO 
1000       %  xdmO 
0       %  ymO 
0 ] ;     %  ydmO 


%  Sampling  time 

%  Missile  system 

%  Target  system 

%  Autopilot 

%  Seeker  head 

[x;xd;y;yd] 


xt(: ,1)=[10000  %  xto 

-2000  %  xdtO 

1000  %  ytO 

0]  ;  %  ydtO 


[x;xd;y;yd] 
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R(l)=sqrt(  (xt(l,l)-xm(l,l)  )  A2+(xt  (3  , 1)  -xm(3  , 1)  )"2)  ; 
time(l)=0;  %  Time 

am(l)=0;amx(l)=0;amy (1)=0;     %  Acceleration 
gamadmis (1)=0;  %  Missile  flight  angle 

beta ( : , 1)=[0;0] ;  %  Seeker  head  angle 

%  BEGIN  SIMULATION  LOOP 

Tf=15;  %  Simulation  time 

kmax=Tf/dt+l; 
for  i=l:kmax-l 

v_mis(i)=sqrt(xm(2, i) A2+xm(4,i) ~2) ;     %  Missile  velocity 

gama_mis(i)=atan2 (xm(4, i) ,xm(2, i) ) ;     %  Missile  vel.  angle 

sigma (1, i)=atan2 ( (xt(3, i)-xm(3, i) ) , (xt(l, i) -xm(l/i) ) ) ; 

beta ( : , i+1) =phis*beta ( : , i) +dels*sigma (i) ; 
%  Input  For  the  Autopilot=betadot*Nr 

gamadmis ( i+1) =phia*gamadmis (i) +dela*Nr*beta (2 , i+1) ; 
%  Missile  Acceleration 

am(i+l)=v_mis (i) *gamadmis (i+1) ; 

if  abs(am(i+l))  >  maxac     %  Check  for  max.  acceleration 
am(i+l) =maxac*sign (am(i+l) ) ; 

end 
%  Acceleration  is  perpendicular  to  gama 

u ( : , i+1) =[-sin(gama_mis (i) ) *am(i+l) ; 
cos (gama_mis (i) ) *am(i+l) ] ; 
%  Update  missile  &  target  states 

xm ( : , i+1) =phim*xm( :  ,  i) +delm*u ( : , i+1) ; 

xt(: ,i+l)=phit*xt(: , i) +delm* [0; 0] ; 

R(i+l)=sqrt(  (xt  (1,  i+1)  -xm(l,  i+1)  )  A2+(xt  (3  ,  i+1)  -xmp,  i+1) )  A2)  ; 
time (i+1) =time (i) +dt; 


%  Check  for  the  min.  distance 

if  abs(R(i+l))  >  abs(R(i)) 
xms=[xm(3/ i-l: i+1) ;xm(4, i-1: i+1) ;xm(l, i-1: i+1) ;xm(2/ i-1: i+1) ] ; 
xts=[xt(3,i-l:i+l) ;xt (4, i-1: i+1) ;xt (1, i-1: i+1) ;xt {2, i-1: i+1) ] ; 

rl  =  interp(xms(: ,1:2) ,xts(: ,1:2) )  ;      %  [Ref.  5] 

r2  =  interp(xms( : ,2:3) ,xts( : ,2:3) ) ; 

rmin=min([rl  r2]); 

break  %  Terminate  if  distance  is  min. 

end 
end 
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%  PLOT  RESULTS 

%  Flight  paths  of  missile  and  target 

plot(xt(l,l:i+l) ,xt(3,l:i+l) ,xm(l, 1:1+1) ,  xm(3 ,1:1+1) ) 

title('Flight  Paths'); 

gtext(['rmin  =  ' /num2str (rmin) ] ) ; pause 

%  Acceleration  Components 

plot (time (1: i) ,u(l,l:i) ,time(l: i) ,u(2, 1: i) ) ; 
title ( 'Acceleration  Components') ;pause; 

%  Missile  acceleration 

plot (time (1: i+1) ,am(l: i+1) ) ; 

title ( 'Missile  Acceleration') ;pause; 

%  Line-of-sight  rate 

plot (time (1: i) , beta (2 , 1:1) ) ; 
title ( 'Line-of-Sight  Rate') ;pause 

%  Line-of-Sight  angle 

plot (time (1 : i) , sigma (1:1)) ; 
title('Line  of  sight' ) ;pause 

%  Range 

plot ( time (1: i+1) ,R(l:i+l) ) 
title ('Range') 
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2.  YYDYDD.M 

%  Missile/Target  engagement  simulation 

%  Control  by  Third  Order  Minimum  Time  Controller 

%  Control  States:  Vertical  Distance  and  derivatives 

%  Missile  acceleration  limited  to  640  ft/sec2 

%  Program  terminates  when  the  range  starts  increasing 

%  Miss  Distance  is  found  by  interpolating  minimum  distance 


clear; 

:clg 

alp=3; 

i 

N=l; 

maxac= 

=640; 

Am=[0 

1  0 

o; 

0 

0  0 

o; 

0 

0  0 

i; 

0 

0  0 

o]; 

Bm=[0 

0;l 

0;0 

o; 

0 

i]; 

Aa=[0 

i; 

0 

-air 

>]; 

%  1/alp  =  Autopilot  time  constant 


%  Missile  matrices 


%  Autopilot  states 


%  Sampling  time 


xm( : ,1:2)=[ 

0 

10; 

%  xmO 

1000 

1000; 

%  xdmO 

. 

0 

o; 

%  ymO 

0 

o]; 

%  ydmO 

xt  ( : ,  1 : 2 )  =  [ 

10000 

9980; 

%  xtO 

-2000 

-2000; 

%  xdtO 

1000 

1000; 

%  ytO 

0 

o]; 

%  ydto 

R(l)=sqrt( (xt(l,l)-xm(l,l) ) ~2+(xt (3 , 1) -xm(3 , 1) ) "2) ; 
R(2)=sqrt( (xt (1 , 2) -xm(l, 2) ) "2+ (xt (3 , 2) -xm(3 , 2) )*2) ; 
time(l:2)=[0  0]; 
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gmis (1)=0; 

gamamis( : , 1:2)=[0  0;0  0]; 

amy ( 1 : 2 ) = [ 0  0 ] ; 

v_mis ( 1 ) =1000 ; aax ( 1 ) =0 ; aay ( 1) =0 ; 

%  SWITCHING  LAW 

ax=0 ;  ay=0 ; 

Vx=(xt(2,l)-xm(2,l) ) ;    Vy=xt(4/l)-xm(4/l) ; 

x=(xt(l,l)-xm(l,l) ) ;    y=xt(3,l)-xm(3,l) ; 

xsl(l)=y;   xs2(l)=Vy;   xs3(l)=0; 

u ( 1) =N*swth_law (xsl , xs2 , xs3 , alp, N) ; 

%BEGIN  SIMULATION  LOOP 

Tf=15; 

kmax=Tf/dt+l; 

for  i=2:kmax-l 

%  Missile  velocity 

v_mis (i)=sqrt (xm(2, i) ~2+xm(4 , i) A2) ; 

gmis (i) =atan2 (xm(4,i),xm(2,i)); 

xsl (i)=(xt (3, i) -xm(3/ i) ) / (v_mis(i) * cos (gmis (i) ) ) ; 

xs2 (i)  =  (xt  (4  ,  i) -xm(4 ,  i)  )  / (v_mis (i) *cos (gmis (i) ) ) ; 

xs3 (i) =-gamamis (2 , i) ;%-amy (i) ;% 

u(i)=N*swth_law(xsl(i) ,xs2 (i) ,xs3 (i) ,alp,N) ; 

%   Input  for  the  autopilot  =  u 

gamamis ( : , i+l) =phia*gamamis ( : , i) +dela*u(i) ; 

%  Missr'.le  acceleration 

am (i+l) =v_mis (i) *gamamis (2 , i+l) ; 

if  abs(am(i+l))  >  maxac 

am (i+l) =maxac*sign(am(i+l) ) ; 

end 

%   Acceleration  components 
amx(i+l) =-am(i+l) * sin (gmis (i) ) ; 
amy(i+l)=  am(i+l) *cos (gmis (i) ) ; 

%   Update  missile  and  target  states 

xm( : , i+l)=phim*xm( : , i) +delm* [  amx(i+l);  amy(i+l)]; 

xt(: /i+l)=phit*xt(: , i) +delm* [0;0] ; 

R(i+l)=sqrt(  (xt  (1,  i+l)  -xm(l,  i+l) )  ~2+(xt(3,  i+l)  -xm(3  ,  i+l) )  ~2)  ; 

time ( i+l) =time(i) +dt; 

if  (abs(R(i+l))  >abs(R(i))) 
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xms=[xm(3,i-l:i+l) ;xm(4, i-1: i+1) ;xm(l, i-l:i+i) ;xm(2, i-i:i+l) ] ; 

xts=[xt(3,i-l:i+l) ;xt(4,i-l:i+l) ;xt (1,1-1:1+1) ;xt(2, i-1: i+1) ] ; 

rl  =  interp(xms( : ,1:2) ,  xts( : ,1:2) ) ;      %  [Ref.  5] 

r2  =  interp(xms( : ,2:3) ,xts( : ,2:3) ) ; 

rmin=min( [rl  r2]); 

break 
end 
end 

%  PLOT  RESULTS 
%  Flight  paths 
plot  (xm(l,  1:3:  i+1)  ,xm(3  , 1:  3  :  i+1)  ,xt(l,  1:3:  i+1)  ,xt  (3, 1:3  :  i+1) ) 

title ( ' Flightpaths ' ) ; 

gtext ( [ 'rmin=' , num2str (rmin) ] ) ;pause; 

%   Control  input 

axis([0  time(i)+0.1  -N-N/10  N+N/10]) 

plot(time(l:2:i) ,u(l:2:i) ,time(l:2:i) ,u(l: 2 : i) , '*' ) ; 
title ('U') ; pause; axis; 

%   Missile  acceleration 

axis([0  time (i) +0.1  -maxac-maxac/25  maxac+maxac/25] ) 

plot (time (1:2 : i+1) ,am(l:2:i+l) ) ; 

title ('Missile  acceleration' ) ; pause; axis; 

%   Missile  velocity  angle 

plot (time (1:2: i) ,gmis (l, 1:2: i) ) ; 

title ( 'Gama' ) ; 

pause; 

%   Missile  velocity 

plot (time (1:2 : i) , v_mis(l:2: i) , time (1:2: i) ,v_mis(l:2 :i) , '*') ; 

title( 'Missile  velocity ') ;pause 

%   y  states 

plot (time (1:2 : i) ,xt(3, 1:2:1) -xm (3, 1:2:1) ) ; 

title (' (Y_target  -  Y_missile  vs  Time  '); 

pause; 
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3.  SIGMDER.M 

%  Missile/Target  engagement  simulation 

%  Control  by  Third  Order  Minimum  Time  Controller 

%  Control  States:  Derivatives  of  Sigma 

%  Missile  acceleration  limited  to  640  ft/sec2 

%  Program  terminates  when  the  range  starts  increasing 

%  Miss  Distance  is  found  by  interpolating  minimum  distance 

clear ;clg 

alp=3; 

N=l; 

maxac=64  0; 


Am=[0  10  0; 
0  0  0  0; 
0  0  0  1; 

0  0  0  0  ]  ; 
Bm=[0  0; 

1  0; 
0  0; 
0  1]; 

Aa= [0    1 ; 

0  -alp] ; 
Ba=[0;alp] ; 
dt=0.01; 
Tf=25; 

kmax=Tf/dt+l; 
%    INITIALIZE  STATES 
xm  ( :  ,  1 )  =  [   0  ; 
1000; 

o; 
o]; 

xt(:,l)=[  10000; 

-2000; 

1000; 

0]; 

xm( : ,2)=[  10; 
1000; 
o; 
o]; 

xt(:,2)=[   9980; 
-2000; 


%  Missile  matrices 


%  Autopilot  states 


%  xmO 

%  xdmO 

%  ymO 

%  ydmO 

%  xtO 

%  xdtO 

%  ytO 

%  ydtO 

%  xmO 

%  xdmO 

%  ymO 

%  ydmO 

%  xtO 

%  xdtO 
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1000;       %  ytO 
0] ;     %  ydto 
R(l)=sqrt((xt(l,l)-xm(l,l)) *2+(xt (3 , 1) -xm(3 , 1) ) A2) ; 
R(2)=sqrt((xt(l,2)-xm(l,2) ) *2+(xt (3 ,2) -xm(3 , 2) )  A2)  ; 
gamamis( : , 1: 2)=[0  0;0  0]; 

am ( 1 : 2 ) = [ 0  0 ] ; amx ( 1 : 2 ) = [ 0  0 ] ;  amy ( 1 : 2 ) = [ 0  0  ] ; 
time(l:2)=[0  dt] ; 
[phim,delm]=c2d(Am,Bm,dt) ; 
[phit]=[phim] ; 
[delt]=[delm] ; 
[phia,dela]=c2d(Aa,Ba,dt) ; 
v_mis ( 1 ) =1 00  0 ; aax ( 1 ) =0 ; aay ( 1 ) =0 ; 
ax=0;  ay=0; 

Vx=(xt(2,l)-xm(2,l)) ;    Vy=xt(4,l)-xm(4,l) ; 
x=(xt(l,l)-xm(l,l) ) ;    y=xt(3,l)-xm(3,l) ; 
sigma (1, l)=atan2 (y,x) ; 
sigma(2, l)=(Vy*x-Vx*y)  /  (xA2+yA2)  ; 

sigma(3, 1) =(ay*x-ax*y) /R(l) "2-2* (x*Vx+y*Vy) * (x*Vy-y*Vx) /R(l) 
M; 

sigma (3 , 1) = (ay*x-ax*y) / (xA2+yA2) -2* (x*Vx+y*Vy) * (x*Vy-y*Vx) / ( 
x"2+yA2) A2; 

xsl (l)=sigma(2/ 1) ;xs2 (l)=sigma(3/ 1) ;xs3 (1)=0; 
u(l)=N*swth_law(xsl,xs2 ,xs3 ,alp,N) ; 
for  i=2:kmax-l 

%  Missile  velocity 

v_mis(i)=sqrt(xm(2, i) A2+xm(4, i) A2) ; 

gama_mis (i) =atan2 (xm(4,i) , xm(2 , i) ) ; 

ay=-amy(i);  ax=-amx(i) ; 

Vy=xt(4,i)-xm(4,i) ;    Vx=(xt (2 , i) -xm(2 , i) ) ; 

y=xt ( 3  ,  i ) -xm ( 3 , i ) ;     x= ( xt ( 1 , i ) -xm (1,1)); 

aay (i)=(-amy (i-1) +v_mis (i-1) *cos (gamamis (1, i) ) *u(i-l) ) ; 

aax(i)=(-amx(i-l)+v_mis(i-l) *s in (gamamis (1, i) ) *u (i-1) ) ; 

sigma (1, i) =atan2 (y,x) ; 

sigma(2, i)=(Vy*x-Vx*y) / (x*2+y~2) ; 

sigma ( 3 ,  i) = (ay*x-ax*y) / (xA2+y~2) -2* (x*Vx+y*Vy) * (x*Vy-y*Vx) / ( 
xA2+yA2) A2; 

xsl (i) =sigma(2, i) ;xs2 (i)=sigma(3/ i) ; 

xs3 (i)=(Vx*ay+x*aay (i) -Vy*ax-y*aax(i) ) /R(i) "2-. . . 

2*(x*Vx+y*Vy)*(x*ay-y*ax)/R(ip4-.  .  . 

2*( (VxA2+x*ax+VyA2+y*ay) * (x*Vy-y*Vx) +(x*ay-y*ax) *(x*Vx+y*Vy) 

)/R(i)M+... 
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8* (x*Vx+y*Vy) A2* (x*Vy-y*Vx) /R(i) A6; 
u(i)=N*swth_law(xsl(i) ,  xs2 (i) ,  xs3 (i) , alp,N) ; 
%   END  SWITCHING  LAW 
%   Input  for  the  autopilot  =  u 
gamamis ( : , i+1) =phia*gamamis ( : , i) +dela*u(i) ; 
%  Missile  acceleration 
am(i+l)=v_mis(i) *gamamis(2, i+1) ; 
if  abs(am(i+l))  >  maxac 

am(i+l)=maxac*sign(am(i+l) ) ; 
end 
%  Acceleration  components 
amx(i+l)=-am(i+l) *sin (gamamis (1, i+1) ) ; 
amy (i+1) =  am (i+1) *cos (gamamis (1, i+1) ) ; 
%   Update  missile  and  target  states 
xm( : , i+l)=phim*xm( : , i) +delm* [  amx(i+l) ;  amy(i+l)]; 
xt ( : , i+1 ) =phit*xt ( :  ,  i ) +delm* [ 0 ; 0 ] ; 

R(i+l)=sqrt(  (xt  (1,  i+1)  -xm(l,  i+1)  )  ~2+(xt  (3  ,  i+1)  -xm(3,  i+1)  )"2)  ; 
time (i+1) =time (i) +dt; 

if  (abs(R(i+l))  >  abs(R(i))) 

xms=[xm(3, i-1: i+1) ;xm(4, i-1: i+1) ;xm(l, i-1: i+1) ;xm(2, i-l:i+l) ]  ; 

xts=[xt(3,i-l:i+l) ;xt (4, i-1: i+1) ;>t (1, i-1: i+1) ;xt (2, i-1: i+1) ] ; 
rl  =  interp(xms( : , 1:2) ,xts( : , 1:2) ) ;    %  [Ref.  5] 
r2  =  interp(xms ( : , 2 : 3) ,xts ( : ,2 : 3) ) ; 
rmin=min([rl  r2 ] ) ; 
break 
end 
end 

%  PLOT  RESULTS 
%     Flight  paths 

plot (xt (1,1: 3: i+1) ,xt(3, 1:3: i+1) , '*' ,xm(l, 1: 3 : i+1) ,xm(3,l:3: 
i+1)) 

title ('Flight  paths'); 
gtext ( [ 'rmin=/ , num2str (rmin) ] ) ; 
pause; 
%   Control  input 

plot (time (1: 2 :i) ,u (1:2:1) ,time(l:2:i)  ,u(l:2  :  i)  ,  '*' )  ; 

title('U') ; 

pause; 
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%   Missile  acceleration 

plot(time(l:2:i+l)  ,310(1:2:1+1)  ,time(l:2:  i+1)  ,3111(1:2:1+1)  ,'*'); 

title ( 'Missile  acceleration'); 

pause; 

%   Missile  velocity  angle 

plot (time(l:2: i) , gamamis(l, 1:2: i) ,time(l:2: i) ,gama_mis(l, 1:2 

:i),'*'); 

title('Gama') ;xlabel ( 'Time' ) ;ylabel ( 'Radian' ) 

pause; 

%   Missile  velocity 

plot (time (1:2: i) ,v_mis (1:2 : i) ,time(l:2:i) , v_mis(l: 2 : i) ,'*'); 

title( 'Missile  velocity'); 

; pause 

%   y  states 

plot (time (l: 2 :i) ,xt(3,l:2:i) -xm(3,l:2: i) ) ; 

title (' (Y_target  -  Y_missile)  vs  Time' ) ;pause;%meta  sigcd6 

%   range 

plot (time, R) ;title( 'range' ) ;pause;%meta  range 
%    sigma 

plot (time (1: i) , sigma (1, 1: i) ) ;title( 'sigma' ) ; 
pause; 
%    sigmadot 

plot  (time  (1:  i)  ,  sigma  (2  ,  l.:i) )  ; title  ( 'sigmadot' )  ; pause 

%meta  sigmadot 

plot (time, amx, time, amy) ; title ( 'comp. ' ) ; pause ;%meta  comp 
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4.  SWTH_LAW.M 

%  Function  calculates  the  sign  of  the  control 

function [s]=swth_law(xl/x2 ,x3 , alp,N) ; 

w=sign( (alp*xl+x2) +( (alp*x2+x3) *abs (alp*x2+x3) / (2*N) ) ) ; 
f=(alp*xl+x2)+w*((alp*x2+x3)"2/(2*N) ) ; 
z=sqrt(abs(f)/N) ; 

s=sign(w*( (N*exp(alp*z) /alp) -(2*N/alp) ) +( (exp(-alp*z) * . 
exp(-w*alp*(alp*x2+x3)/N) ) * (x3+w*N/alp) ) ) ; 
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